box::use(
  ggdist[...],
  dplyr[select, mutate, across, group_by, summarise, arrange, case_when, filter, if_else, distinct, ungroup, row_number, add_row, slice, n],
  patchwork[...],
  data.table[setDT, rbindlist],
  modelsummary[...],
  forcats[as_factor, fct_relevel],
  estimatr[...],
  distributional[...],
  gghalves[...],
  stringr[...],
  qte[...],
  tidyr[pivot_wider, separate, pivot_longer],
  ./helpers/theme_helpers[conjoint_reshape, plexi, format_si, percent, clean_cj_amce],
  rlang[...],
  ggplot2[...],
  ./helpers/descriptives,
  pdftools[...]
)

# ------------------------------- cleaned data ------------------------------- #

df <- readRDS("./data/df_ordering.RDS")


# ---------------------------------------------------------------------------- #
#                                Figures 1 and 2; A4                           #
# ---------------------------------------------------------------------------- #

# ---------------------------------------------------------------------------- #
#                                     Fig 2                                    #
# ---------------------------------------------------------------------------- #

fit_models_cog_burden_imp_first <- function(df, dvs) {

  models <- vector("list", length(dvs))

  for (x in 1:length(dvs)) {
    
    models[[x]] <- tidy(lm_robust(df, formula = as.formula(paste(dvs[x],"~ econ_moral_first")))) %>% 
    add_row(term = "Baseline: Important attributes \nnot at top", estimate = 0, outcome = dvs[x], .before = c(2, 4, 6))
  }

  models <- rbindlist(models, idcol=T)

  models <- models[term != "(Intercept)"]

  models <- models %>% 
    mutate(term = ifelse(term == "econ_moral_firstImportant attributes first", 
                       gsub("econ_moral_first", "", term), term))
}


cog_burden_imp_first <- df %>% 
  select(uuid, rand_type = CJ_ATT_ORDERING, Q16r1, Q16r2, Q16r4, BLOCKORDERc1, BLOCKORDERc2, BLOCKORDERc3, BLOCKORDERc4, BLOCKORDERc5) %>% 
  mutate(rand_type = as_factor(if_else(rand_type == 1, "Full Randomization", "Block Randomization"))) %>% 
  mutate(across(Q16r1:Q16r4, ~ str_remove(., "\\("))) %>% 
  mutate(across(Q16r1:Q16r4, ~ str_remove(., "\\)"))) %>%
  mutate(across(Q16r1:Q16r4, ~ as.numeric(.))) %>% 
  filter(rand_type == "Block Randomization") %>% 
  mutate(econ_moral_first = case_when(
    BLOCKORDERc1 == 1 | BLOCKORDERc3 == 1 ~ "Important attributes first",
    TRUE ~ "Else"
  ))

 
dvs <- c("Q16r1", "Q16r2", "Q16r4")

models <- fit_models_cog_burden_imp_first(cog_burden_imp_first, dvs)  

models_q16 <- models %>% 
mutate(outcome = case_when(outcome == "Q16r1" ~ "I often did not know what to answer.", 
                           outcome == "Q16r2" ~ "It was easy to make a decision between the requests for arms exports.",
                           outcome == "Q16r4" ~ "I was missing important information in the tables to make a proper decision."
                           )) %>% 
                           slice(1:6)


models_q16 <- models_q16 %>% 
 mutate(order = row_number()*-1) 

models_q16$order <- factor(models_q16$order) 


p <- ggplot(models_q16, aes(y=order,x=estimate)) + xlim(-0.5,0.5)  + geom_point(size=2)+
  geom_segment(data=models_q16, aes(x=conf.low,xend=conf.high, y=order,yend=order)) +
  geom_text(aes(label=ifelse(estimate != 0, round(estimate, digits = 3),"")), vjust = 0, nudge_y = 0.3, size = 5) +
  facet_wrap(outcome~., scales="free_y", ncol=1, labeller = label_wrap_gen(multi_line = TRUE)) +
  scale_y_discrete(breaks=models_q16$order,labels=models_q16$term)+ ylab("") +
  theme(strip.text.y=element_text(angle=0))  +
  xlab("\nATE") +
  geom_vline(xintercept = 0,size=.5,colour="black",linetype="dotted") + ggtitle("How much do you agree or disagree with the following statements?") + plexi() 

ggsave(p, filename = "./figures/fig_2_cog_burden_imp_first_q16_ATEs.pdf", width = 10, height = 6, device = cairo_pdf)

pdf_convert("./figures/fig_2_cog_burden_imp_first_q16_ATEs.pdf", filenames = "./figures/fig_2_cog_burden_imp_first_q16_ATEs.tiff", dpi = 300, format = "tiff")

# ---------------------------------------------------------------------------- #
#                                     Fig 1                                    #
# ---------------------------------------------------------------------------- #


cog_burden_imp_first_timing <- df %>% 
  select(uuid, rand_type = CJ_ATT_ORDERING, CJ_1_timer, BLOCKORDERc1, BLOCKORDERc2, BLOCKORDERc3, BLOCKORDERc4, BLOCKORDERc5) %>% 
  filter(CJ_1_timer <= 600) %>% 
  mutate(rand_type = as_factor(if_else(rand_type == 1, "Full Randomization", "Block Randomization"))) %>% 
  filter(rand_type == "Block Randomization") %>% 
  mutate(econ_moral_first = case_when(
    BLOCKORDERc1 == 1 | BLOCKORDERc3 == 1 ~ "Important attributes first",
    TRUE ~ "Else"
  )) %>% 
  mutate(binary_treat = if_else(econ_moral_first == "Else", 0, 1)) %>% 
  group_by(econ_moral_first) %>% 
  mutate(n = n()) %>% 
  ungroup()


models_timing <- fit_models_cog_burden_imp_first(cog_burden_imp_first_timing, dvs = "CJ_1_timer")

p1 <- ggplot() +
  stat_pointinterval(aes(x = models_timing$term[1], y = models_timing$estimate[1])) +
  stat_halfeye(
    aes(
      x = models_timing$term[2], ydist = dist_student_t(df = models_timing$df[2], mu = models_timing$estimate[2], sigma = models_timing$std.error),
      fill = stat(cut_cdf_qi(cdf, .width = c(.66, .95, 1)))
    ),
    alpha = 0.7, slab_colour = "#0f7ba6"
  ) +
  scale_fill_manual(values = c("#0f7ba6", "skyblue", "gray85")) +
  guides(fill = "none") +
  labs(x = "", y = "ATE \n(in seconds)") +
  plexi() +
  scale_y_continuous(position = "right", breaks = seq(-6, 6, by = 1)) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  geom_hline(yintercept = models_timing$estimate[2], linetype = "dashed") +
  theme(panel.grid.major.y = element_line(colour = "grey92"), 
        panel.grid.major.x = element_blank())+
  theme(axis.title.y.right = element_text(angle = 0, hjust = 0), axis.text.y = element_text(size = 14), axis.text.x = element_text(size = 14), axis.title=element_text(size=14))

p1

p2 <- cog_burden_imp_first_timing %>% 
  mutate(econ_moral_first = if_else(econ_moral_first == "Else", paste0("Important attributes \n not at top (N = ", distinct(., n)$n[2], ")"), paste0("Important attributes first \n (N = ", distinct(., n)$n[1], ")"))) %>% 
  ggplot(aes(x=econ_moral_first, y=CJ_1_timer)) +
  stat_halfeye(
    adjust = 1, 
    width = .6, 
    .width = 0, 
    justification = -.2, 
    point_colour = NA,
    slab_colour = "#0f7ba6",
    slab_fill = "skyblue",
    slab_alpha = 0.7
  ) + 
  geom_boxplot(
    width = .15, 
    outlier.shape = NA
  ) +
  gghalves::geom_half_point(
    side = "l", 
    range_scale = .5, 
    alpha = .1
  ) +
  coord_cartesian(xlim = c(1.2, NA), clip = "off") +
  labs(y = "Completion Time \nTask 1 (in seconds)", x = "") +
  scale_y_continuous(breaks = seq(0, 600, by = 50)) + plexi() +
  theme(panel.grid.major.y = element_line(colour = "grey92"), 
        panel.grid.major.x = element_blank()) +
  theme(axis.title.y.left = element_text(angle = 0, hjust = 1), axis.text.x = element_text(size = 14), axis.text.y = element_text(size = 14), axis.title=element_text(size=14))
  
p2

p2 + p1

ggsave(filename = "./figures/fig_1_cog_burden_imp_first_times.pdf", width = 17, height = 6, device = cairo_pdf)

pdf_convert("./figures/fig_1_cog_burden_imp_first_times.pdf", filenames = "./figures/fig_1_cog_burden_imp_first_times.tiff", dpi = 300, format = "tiff")

# ---------------------------------------------------------------------------- #
#                                     Fig A1                                   #
# ---------------------------------------------------------------------------- #

qtes <- ci.qte(CJ_1_timer ~ binary_treat, x = NULL,
   data=cog_burden_imp_first_timing,
   probs=seq(0.05,0.95,0.05), se=T, iters=10)

p = ggqte(qtes, xbreaks = seq(0,1,0.1), ybreaks = seq(-20, 20, 2), ylab = "QTE \n(time in seconds)") + plexi() +
  theme(axis.title.y.left = element_text(angle = 0, hjust = 1), axis.text.x = element_text(size = 14), axis.text.y = element_text(size = 14), axis.title=element_text(size=14))
  
ggsave(p, filename = "./figures/fig_a1_cog_burden_imp_first_time_qte.pdf", width = 7, height = 6, device = cairo_pdf)

pdf_convert("./figures/fig_a1_cog_burden_imp_first_time_qte.pdf", filenames = "./figures/fig_a1_cog_burden_imp_first_time_qte.tiff", dpi = 300, format = "tiff")

# ---------------------------------------------------------------------------- #
#                                Figures 3 and 4; A5                           #
# ---------------------------------------------------------------------------- #

# ---------------------------------------------------------------------------- #
#                                     Fig 4                                    #
# ---------------------------------------------------------------------------- #

fit_models <- function(df, dvs) {

  models <- vector("list", length(dvs))

  for (x in 1:length(dvs)) {
    
    models[[x]] <- tidy(lm_robust(df, formula = as.formula(paste(dvs[x],"~ rand_type")))) %>% 
    add_row(term = "Baseline: Full Randomization", estimate = 0, outcome = dvs[x], .before = c(2, 4, 6))
  }

  models <- rbindlist(models, idcol=T)

  models <- models[term != "(Intercept)"]

  models <- models %>% 
    mutate(term = ifelse(term == "rand_typeBlock Randomization", 
                       gsub("rand_type", "", term), term))
}


cog_burden_block <- df %>% 
  select(uuid, rand_type = CJ_ATT_ORDERING, Q16r1, Q16r2, Q16r4) %>% 
  mutate(rand_type = as_factor(if_else(rand_type == 1, "Full Randomization", "Block Randomization"))) %>% 
  mutate(across(Q16r1:Q16r4, ~ str_remove(., "\\("))) %>% 
  mutate(across(Q16r1:Q16r4, ~ str_remove(., "\\)"))) %>%
  mutate(across(Q16r1:Q16r4, ~ as.numeric(.)))

  
dvs <- c("Q16r1", "Q16r2", "Q16r4")

models <- fit_models(cog_burden_block, dvs)

models_q16 <- models %>% 
mutate(outcome = case_when(outcome == "Q16r1" ~ "I often did not know what to answer.", 
                           outcome == "Q16r2" ~ "It was easy to make a decision between the requests for arms exports.",
                           outcome == "Q16r4" ~ "I was missing important information in the tables to make a proper decision."
                           )) %>% 
                           slice(1:6)


models_q16 <- models_q16 %>% 
 mutate(order = row_number()*-1) 

models_q16$order <- factor(models_q16$order) 

p <- ggplot(models_q16, aes(y=order,x=estimate)) + xlim(-0.5,0.5) + geom_point(size=2)+ plexi() +
  geom_segment(data=models_q16, aes(x=conf.low,xend=conf.high, y=order,yend=order)) +
  geom_text(aes(label=ifelse(estimate != 0, round(estimate, digits = 3),"")), vjust = 0, nudge_y = 0.3, size = 5) +
  facet_wrap(outcome~., scales="free_y", ncol=1) +
  scale_y_discrete(breaks=models_q16$order,labels=models_q16$term)+ ylab("") +
  theme(strip.text.y=element_text(angle=0))  +
  xlab("\nATE") + theme(
    text = element_text(family = "CMU Serif", size = 16, colour="black")) +
  geom_vline(xintercept = 0,size=.5,colour="black",linetype="dotted") + ggtitle("How much do you agree or disagree with the following statements?")

ggsave(p, filename = "./figures/fig_4_cog_burden_block_q16_ATEs.pdf", width = 11, height = 5, device = cairo_pdf)

pdf_convert("./figures/fig_4_cog_burden_block_q16_ATEs.pdf", filenames = "./figures/fig_4_cog_burden_block_q16_ATEs.tiff", dpi = 300, format = "tiff")

# ---------------------------------------------------------------------------- #
#                                     Fig 3                                    #
# ---------------------------------------------------------------------------- #

cog_burden_block_timing <- df %>% 
  select(uuid, rand_type = CJ_ATT_ORDERING, CJ_1_timer) %>% 
  filter(CJ_1_timer <= 600) %>% 
  mutate(rand_type = as_factor(if_else(rand_type == 1, "Full Randomization", "Block Randomization"))) %>% 
  mutate(binary_treat = if_else(rand_type == "Full Randomization", 0, 1)) %>% 
  group_by(rand_type) %>%
  mutate(n = n()) %>% 
  ungroup()


models_timing <- fit_models(cog_burden_block_timing, dvs = "CJ_1_timer")

p1 <- ggplot() +
  stat_pointinterval(aes(x = models_timing$term[1], y = models_timing$estimate[1])) +
  stat_halfeye(
    aes(
      x = models_timing$term[2], ydist = dist_student_t(df = models_timing$df[2], mu = models_timing$estimate[2], sigma = models_timing$std.error),
      fill = stat(cut_cdf_qi(cdf, .width = c(.66, .95, 1)))
    ),
    alpha = 0.7, slab_colour = "#0f7ba6"
  ) +
  scale_fill_manual(values = c("#0f7ba6", "skyblue", "gray85")) +
  guides(fill = "none") +
  labs(x = "", y = "ATE \n(in seconds)") +
  plexi() +
  scale_y_continuous(position = "right", breaks = seq(-4, 6, by = 0.5)) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  geom_hline(yintercept = models_timing$estimate[2], linetype = "dashed") +
  theme(panel.grid.major.y = element_line(colour = "grey92"), 
        panel.grid.major.x = element_blank())+
  theme(axis.title.y.right = element_text(angle = 0, hjust = 0),axis.text.x = element_text(size = 14), axis.text.y = element_text(size = 14), axis.title=element_text(size=14))


p2 <- cog_burden_block_timing %>% 
  mutate(rand_type = if_else(rand_type == "Full Randomization", paste0("Full Randomization \n (N = ", distinct(., n)$n[1], ")"), paste0("Block Randomization \n (N = ", distinct(., n)$n[2], ")"))) %>% 
  ggplot(aes(x=rand_type, y=CJ_1_timer)) +
  stat_halfeye(
    adjust = 1, 
    width = .6, 
    .width = 0, 
    justification = -.2, 
    point_colour = NA,
    slab_colour = "#0f7ba6",
    slab_fill = "skyblue",
    slab_alpha = 0.7
  ) + 
  geom_boxplot(
    width = .15, 
    outlier.shape = NA
  ) +
  gghalves::geom_half_point(
    side = "l", 
    range_scale = .5, 
    alpha = .1
  ) +
  coord_cartesian(xlim = c(1.2, NA), clip = "off") +
  labs(y = "Completion Time \nTask 1 (in seconds)", x = "") +
  scale_y_continuous(breaks = seq(0, 600, by = 50)) + plexi() +
  theme(panel.grid.major.y = element_line(colour = "grey92"), 
        panel.grid.major.x = element_blank()) +
  theme(axis.title.y.left = element_text(angle = 0, hjust = 1), axis.text.x = element_text(size = 14), axis.text.y = element_text(size = 14), axis.title=element_text(size=14))
  
p2 + p1

ggsave(filename = "./figures/fig_3_cog_burden_block_times.pdf", width = 17, height = 6, device = cairo_pdf)

pdf_convert("./figures/fig_3_cog_burden_block_times.pdf", filenames = "./figures/fig_3_cog_burden_block_times.tiff", dpi = 300, format = "tiff")

# ---------------------------------------------------------------------------- #
#                                     Fig A2                                   #
# ---------------------------------------------------------------------------- #

qtes <- ci.qte(CJ_1_timer ~ binary_treat, x = NULL,
   data=cog_burden_block_timing,
   probs=seq(0.05,0.95,0.05), se=T, iters=10)

p = ggqte(qtes, xbreaks = seq(0,1,0.1), ybreaks = seq(-20, 25, 2), ylab = "QTE \n(time in seconds)") + plexi() +
  theme(axis.title.y.left = element_text(angle = 0, hjust = 1), axis.text.x = element_text(size = 14), axis.text.y = element_text(size = 14), axis.title=element_text(size=14))
  
ggsave(p, filename = "./figures/fig_a2_cog_burden_block_time_qte.pdf", width = 7, height = 6, device = cairo_pdf)

pdf_convert("./figures/fig_a2_cog_burden_block_time_qte.pdf", filenames = "./figures/fig_a2_cog_burden_block_time_qte.tiff", dpi = 300, format = "tiff")

# ---------------------------------------------------------------------------- #
#                                   Figure 6                                   #
# ---------------------------------------------------------------------------- #

df_conjoint <- conjoint_reshape(df)

att_names <- c(
  `ATT1` = "Regime Type",
  `ATT2` = "Human Rights \nSituation",
  `ATT3` = "Conflict Situation",
  `ATT4` = "Security Partner",
  `ATT5` = "Trade Relationship",
  `ATT6` = "Economic Value",
  `ATT7` = "Jobs at Risk",
  `ATT8` = "Weapon Type",
  `ATT9` = "Weapon Exports \nof other Countries"
)

df_h4 <- df_conjoint %>% 
  select(uuid, CH, CJ,
  y, 
  ATT1:ATT9,
  CONJOINTORDERc1:CONJOINTORDERc9)


ordering_effects <- function(df, treatments) {
    
   out <- vector(mode = "list", length = length(treatments))

   for (i in seq_along(treatments)) {

    h4 <- df %>% 
      filter(!!sym(paste0("CONJOINTORDERc", i)) == 1 | !!sym(paste0("CONJOINTORDERc", i)) %in% 2:8) %>% 
      mutate(first_rest = if_else(!!sym(paste0("CONJOINTORDERc", i)) == 1, "first", "rest")) 

    first <- tidy(lm_robust(data = h4  %>% filter(!!sym(paste0("CONJOINTORDERc", i)) == 1), formula(paste0("y ~ ", treatments[i])), clusters = uuid, se_type = "stata")) %>% 
      filter(term != "(Intercept)") %>% 
      select(term, estimate, conf.low, conf.high)

    rest <- tidy(lm_robust(data = h4  %>% filter(!!sym(paste0("CONJOINTORDERc", i)) %in% 2:8), formula(paste0("y ~ ", treatments[i])), clusters = uuid, se_type = "stata")) %>% 
      filter(term != "(Intercept)") %>% 
      select(term, estimate, conf.low, conf.high)

      int <- tidy(lm_robust(data = h4, formula(paste0("y ~ ", treatments[i], "*first_rest")), clusters = uuid, se_type = "stata")) %>% 
      filter(str_detect(term, ":") == TRUE) %>%  
      mutate(term = str_remove(term, ":")) %>% 
      mutate(term = str_remove(term, "first_lastlast"))  %>% 
      mutate(attribute = str_extract(term, "ATT[1-9]")) %>% 
      mutate(term = str_remove(term, "ATT[1-9]"))

    factors <- pivot_longer(h4, treatments[i], names_to="attribute", values_to="term") %>% 
      select(term, attribute) %>% 
      distinct() %>% 
      arrange(attribute, term)

    int <- merge(int, factors,  by = c("attribute", "term"), all = TRUE)
    int[["estimate"]][is.na(int[["estimate"]])] <- 0
    int <- setDT(int)
    int[estimate == 0, term := paste("(Baseline =", term, sep = " ")]
    int[estimate == 0, term := paste(term, ")", sep = "")]

    
    first$attribute <- (str_extract(first$term, "ATT[1-9]"))
    first$term <- gsub("ATT[1-9]", "", first$term)

    first <- merge(first, factors, by = c("attribute", "term"), all = TRUE)
    first[["estimate"]][is.na(first[["estimate"]])] <- 0
    first <- setDT(first)
    first[estimate == 0, term := paste("(Baseline =", term, sep = " ")]
    first[estimate == 0, term := paste(term, ")", sep = "")]


    first <- first %>%
    mutate(abs = abs(estimate)) %>%
    group_by(attribute) %>%
    arrange(abs, .by_group = TRUE) %>%
    ungroup() %>%
    mutate(order = row_number() * -1) %>%
    select(-abs)

    rest$attribute <- (str_extract(rest$term, "ATT[1-9]"))
    rest$term <- gsub("ATT[1-9]", "", rest$term)

    rest <- merge(rest, factors, by = c("attribute", "term"), all = TRUE)
    rest[["estimate"]][is.na(rest[["estimate"]])] <- 0
    rest <- setDT(rest)
    rest[estimate == 0, term := paste("(Baseline =", term, sep = " ")]
    rest[estimate == 0, term := paste(term, ")", sep = "")]


    rest <- rest %>%
    arrange(match(rest$term, first$term)) %>%
    mutate(order = row_number() * -1)

    int <- int %>%
    arrange(match(int$term, first$term)) %>%
    mutate(order = row_number() * -1)

    first <- first %>%
      mutate(cond = "First")

    rest <- rest %>%
      mutate(cond = "Position 2-8")

    m_int_cond <- rbind(first, rest)  
    
    m_int_cond$cond <- factor(m_int_cond$cond)

    m_int_cond$order <- factor(m_int_cond$order)

    int$order <- factor(int$order)

    
p <- ggplot(m_int_cond, aes(y = order, x = estimate, shape = cond)) +
  xlim(-0.25, 0.15) +
  geom_point(size = 2, position = position_dodge(width = -1)) +
  plexi() +
  geom_linerange(data = m_int_cond, aes(xmin = conf.low, xmax = conf.high, y = order), position = position_dodge(width = -1), key_glyph = "path") +
  facet_grid(attribute ~ ., scales = "free_y", space = "free", labeller = as_labeller(att_names)) +
  scale_y_discrete(breaks = m_int_cond$order, labels = m_int_cond$term) +
  ylab("") +
  theme(strip.text.y = element_text(angle = 0, size = 16), axis.text = element_text(size = 16)) +
  xlab("") +
  geom_vline(xintercept = 0, size = .5, colour = "black", linetype = "dotted") +
  theme(
    strip.background = element_blank(),
    strip.text.y = element_blank(),
    plot.margin = margin(0, 0, 0, 0, "pt")
  )


p3 <- ggplot(int, aes(y = order, x = estimate)) +
  xlim(-0.1, 0.075) +
  geom_point(size = 2) +
  plexi() +
  geom_segment(data = int, aes(x = conf.low, xend = conf.high, y = order, yend = order)) +
  facet_grid(attribute ~ ., scales = "free_y", space = "free", labeller = as_labeller(att_names)) +
  theme(strip.text.y = element_text(angle = 0, size = 16), axis.text = element_text(size = 16)) +
  geom_vline(xintercept = 0, size = .5, colour = "black", linetype = "dotted") +
  scale_y_discrete(breaks = int$order, labels = int$term) +
  ylab("") +
  theme(strip.text.y = element_text(angle = 0)) +
  xlab("") +
  theme(
    axis.text.y = element_blank(),
    axis.ticks.y = element_blank(),
    axis.title.y = element_blank(),
    plot.margin = margin(0, 0, 0, 0, "pt")
  ) 

if (i == 1) {
    p_full <- p + ggtitle("Sub-Group AMCE") +
    theme(legend.position="none") + p3 + ggtitle("Interaction")
} else if (i <= 8 & i >= 2) {
     p_full <- p +
    theme(legend.position="none") + p3
  } else {
     p_full <- p + theme(
  legend.title = element_blank(),
  legend.direction = "vertical",
  legend.position = "bottom",
  legend.justification = "left",
  legend.text = element_text(size = 14),
  legend.key.height = unit(1, "cm"),
  legend.spacing.y = unit(-0.8, "cm")
) + p3 

  }

    out[[i]] <- p_full

   }

  return(out)

}


plots <- ordering_effects(df = df_h4, treatments = names(att_names))

ordering_effect_plots <- wrap_plots(plots, ncol = 1)

ggsave(filename = "./figures/fig_6_ordering_effect_plots.pdf", ordering_effect_plots, width = 12, height = 16, device = cairo_pdf)

pdf_convert("./figures/fig_6_ordering_effect_plots.pdf", filenames = "./figures/fig_6_ordering_effect_plots.tiff", dpi = 300, format = "tiff")

# ---------------------------------------------------------------------------- #
#                                   Figure 7                                   #
# ---------------------------------------------------------------------------- #

df_conjoint <- df_conjoint %>% 
  mutate(rand_type = as_factor(if_else(CJ_ATT_ORDERING == 1, "Full Randomization", "Block Randomization"))) %>% 
  mutate(rand_type = fct_relevel(rand_type))

df_cj_full <- df_conjoint %>% 
  filter(rand_type == "Full Randomization")

df_cj_block <- df_conjoint %>% 
  filter(rand_type == "Block Randomization")


m2 <-   tidy(lm_robust(data = df_cj_full, y ~ ATT1 + ATT2 + ATT3+ ATT4 + ATT5 + ATT6 + ATT7 + ATT8 + ATT9, clusters = uuid, se_type = "stata"))

m2 <- clean_cj_amce(df = m2, input = df_conjoint)

m3 <-  tidy(lm_robust(data = df_cj_block, y ~ ATT1 + ATT2 + ATT3+ ATT4 + ATT5 + ATT6 + ATT7 + ATT8 + ATT9, clusters = uuid, se_type = "stata"))

m3 <- clean_cj_amce(df = m3, input = df_conjoint)

m3 <- m3 %>% 
  arrange(match(m2$term, m3$term)) %>% 
  mutate(order = row_number()*-1)

# Differences in CAMCEs

m4 <- tidy(lm_robust(data = df_conjoint, y ~ ATT1*rand_type + ATT2*rand_type + ATT3*rand_type + 
                       ATT4*rand_type + ATT5*rand_type + ATT6*rand_type + ATT7*rand_type + ATT8*rand_type + ATT9*rand_type, clusters = uuid, se_type = "stata"))



# m4 <-  tidy(lm_robust(data = cj_tidy, y ~ ATT1*de + ATT2*de + ATT3*de 
#                       + ATT4*de + ATT5*de + ATT6*de + ATT7*de + ATT8*de + ATT9*de,  cluster = uuid, method = "stata"))

setDT(m4)
m4 <- m4[term != "(Intercept)"]
m4$attribute <- (str_extract(m4$term, "ATT[1-9]"))
m4$term <- gsub("ATT[1-9]", "", m4$term)
m4 <- m4[grepl("rand_typeBlock Randomization:", m4$term) | term == "Not democratically elected:rand_typeBlock Randomization",]
m4$term <- gsub("rand_typeBlock Randomization:", "", m4$term)
m4$term <- gsub("Not democratically elected:rand_typeBlock Randomization", "Not democratically elected", m4$term)


factors <- pivot_longer(df_conjoint, ATT1:ATT9, names_to="attribute", values_to="term") %>% 
  select(term,attribute) %>% 
  distinct() %>% 
  arrange(attribute, term)


m4 <- merge(m4,factors,  by = c("attribute", "term"), all = TRUE)
m4[["estimate"]][is.na(m4[["estimate"]])] <- 0
m4[estimate == 0, term := paste("(Baseline =", term, sep = " ")]
m4[estimate == 0, term := paste(term, ")", sep = "")]


m4 <- m4[match(m2$term, m4$term),]

m4 <- m4 %>% 
  mutate(order = row_number()*-1)

m4$order <- factor(m4$order)



p2 <- ggplot(m2, aes(y=order,x=estimate)) + xlim(-0.25,0.2) + geom_point(size=2)+ plexi() +
  geom_segment(data=m2, aes(x=conf.low,xend=conf.high, y=order,yend=order)) +
  geom_text(aes(label=ifelse(estimate != 0, round(estimate, digits = 3),"")), vjust = 1, nudge_y = 0.3,  hjust =  -.5, size = 5) + 
  facet_grid(attribute~., scales="free_y", space="free", labeller = as_labeller(att_names)) +
  scale_y_discrete(breaks=m2$order,labels=m2$term)+ ylab("") +
  xlab("") +
  geom_vline(xintercept = 0,size=.5,colour="black",linetype="dotted") + theme(
    strip.background = element_blank(),
    strip.text.y = element_blank(),
    axis.text.y = element_text(size = 14, hjust=0.95,vjust=0.2),
    axis.text.x = element_text(size = 18),
    axis.title = element_text(size = 18)
  )



p3 <- ggplot(m3, aes(y=order,x=estimate)) + xlim(-0.25,0.2) + geom_point(size=2)+ plexi() +
  geom_segment(data=m3, aes(x=conf.low,xend=conf.high, y=order,yend=order)) +
  geom_text(aes(label=ifelse(estimate != 0, round(estimate, digits = 3),"")), vjust = 1, nudge_y = 0.3,  hjust =  -.5, size = 5) +
  facet_grid(attribute~., scales="free_y", space="free", labeller = as_labeller(att_names)) +
  scale_y_discrete(breaks=m3$order,labels=m3$term)+ ylab("") +
  theme(strip.text.y=element_text(angle=0))  +
  xlab("\nAMCE") +
  geom_vline(xintercept = 0,size=.5,colour="black",linetype="dotted") + theme(axis.text.y = element_blank(),
                                                                              axis.ticks.y = element_blank(),
                                                                              axis.title.y = element_blank() ) +
  theme(
    strip.background = element_blank(),
    strip.text.y = element_blank(),
    axis.title = element_text(size = 18),
    axis.text.x = element_text(size = 18)
  )




p4 <- ggplot(m4, aes(y=order,x=estimate)) + xlim(-0.05,0.06) + geom_point(size=2)+ plexi() +
  geom_segment(data=m4, aes(x=conf.low,xend=conf.high, y=order,yend=order)) +
  geom_text(aes(label=ifelse(estimate != 0, round(estimate, digits = 3),"")), vjust = 1, nudge_y = 0.3, hjust = -1.3, size = 5) +
  facet_grid(attribute~., scales="free_y", space="free", labeller = as_labeller(att_names)) +
  scale_y_discrete(breaks=m4$order,labels=m4$term)+ ylab("") +
  theme(strip.text.y=element_text(angle=0, size = 18))  +
  xlab("") +
  geom_vline(xintercept = 0,size=.5,colour="black",linetype="dotted") + theme(axis.text.y = element_blank(),
                                                                              axis.ticks.y = element_blank(),
                                                                              axis.title.y = element_blank(),
    axis.text.x = element_text(size = 18),
    axis.title = element_text(size = 18))



p2 + ggtitle("Full Randomization \n(N = 2260)") + p3 + ggtitle("Block Randomization \n(N = 4357)") + p4 + ggtitle("Difference in CAMCEs \n(Block-Full)")

ggsave(filename = "./figures/fig_7_amce_full_block.pdf", width = 17, height = 10, device = cairo_pdf)

pdf_convert("./figures/fig_7_amce_full_block.pdf", filenames = "./figures/fig_7_amce_full_block.tiff", dpi = 300, format = "tiff")
# ---------------------------------------------------------------------------- #
#                                 descriptives                                 #
# ---------------------------------------------------------------------------- #

descriptives$make_descriptive_tables()